Discrete transverse superconducting modes in nano-cylinders 
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Spatial variation in the superconducting order parameter becomes significant when the system 
is confined at dimensions well below the typical superconducting coherence length. Motivated by 
recent experimental success in growing single-crystal metallic nanorods, we study quantum confine- 
ment effects on superconductivity in a cylindrical nanowire in the clean limit. For large diameters, 
where the transverse level spacing is smaller than superconducting order parameter, the usual ap- 
proximations of Ginzburg-Landau theory are recovered. However, under external magnetic field 
the order parameter develops a spatial variation much stronger than that predicted by Ginzburg- 
Landau theory, and gapless superconductivity is obtained above a certain field strength. At small 
diameters, the discrete nature of the transverse modes produces significant spatial variations in the 
order parameter with increased average magnitude and multiple shoulders in the magnetic response. 

PACS numbers: 74.78.-w, 85.25.-j 



I. INTRODUCTION 

Recent developments in nanofabrication techniques al- 
low access to new physical regimes where various intrin- 
sic order parameters interact with a tuneable confining 
environment. Such order parameters cover an array of 
diverse physical systems, such as ferromagnetsi, quan- 
tum dots 2 , molecular electronics 3 , photonic crystals^, 
superconductors^Sii etc. Quantum-confined supercon- 
ductivity is particularly interesting for its macroscopic 
quantum nature; its well-understood microscopic mech- 
anism can also serve as a platform for studies of other 
many-body nanoscopic quantum confinement effects. 

Since the advent of BCS theory, a great deal has been 
understood in both the microscopic and phenomenologi- 
cal aspects of superconductivity for conventional phonon- 
mediated pairing systems. Theories have been imme- 
diately applied with great success to small supercon- 
ductors?'- of various geometrjaiiiS. Particularly useful 
has been the phenomenological Ginzburg-Landau theory, 
which describes superconductivity directly in terms of the 
superconducting order parameters without appealing to 
an underlying electronic basis. 

Partly due to its enormous success, however, the 
Ginzburg-Landau theory has sometimes been applied be- 
yond its strict regime of validity, especially in systems of 
small size. Part of the justification for this has been 
that experimental samples have often been disordered 
or polycrystalline, in which case confinement effects are 
less pronounced than they are for single crystals. Well- 
established dirty-limit theories for small superconduc- 
tors with strong disorder.^ describe fascinating physics, 
such as gapless superconductivity, down to nanometer 
scales. In the clean limit of microscopic BCS theory, 
where the mean free path I is longer than the coher- 
ence length £o {(- 3> £o) 5 the superconducting behavior 
of small samples is very different from that in the dirty 
limit (I <C £o)^^- Recent experimental techniques for 
producing high quality single crystalline nanostructures 



through electrodeposition into extended nanoporesi 3 - de- 
mands a re-examination of the phenomenology of super- 
conductivity in such systems, working from microscopic 
theories. We specifically aim to investigate the often- 
overlooked spatial structures of the superconducting or- 
der parameter in the confined direction by directly solv- 
ing the Bogoliubov-de Gennes equation and comparing 
with other theories. 

For the last few decades, work on one-dimensional su- 
perconductivity has mostly focused on fluctuation ef- 
fects^i^. These treatments assume featureless trans- 
verse superconducting modes within superconducting 
nanowires and instead concentrate on the physics of 
phase slips in the axial direction. Here we complement 
these previous approaches by considering the effects of 
transverse quantum confinement on the spatial variation 
of superconducting order parameter, with consequences 
for the quasi-particle excitation spectrum and the mag- 
netic response. 



II. FORMALISM 

We consider a superconducting cylinder with a radius 
R smaller than the penetration depth A, but much larger 
than the atomic scale, so that we can describe the sys- 
tem with a continuum basis. The Bogoliubov-de Gennes 
(BdG) equations are 
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where A is the order parameter and Hq is the Hamilto- 
nian for electrons, 
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Here m* is the band electron mass, fi is the chemical 
potential, \xb is the Bohr magneton, a is the Pauli spin 
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matrix and H is the external magnetic field. Ho and its 
complex conjugate Hq act on the time-reversed electrons 
in Cooper pairs. The state [u, v] represents the ampli- 
tudes of the pair of electrons which interact with each 
other via the pairing interaction parametrized by the su- 
perconducting order parameter A. For an axial magnetic 
field, H = Hz, the vector potential in the Coulomb gauge 
is: 



A = -rHO, 
2 



(3) 



where 9 is an unit vector along the azimuthal direction. 
We assume that the radius of the nanorod is sufficiently 
below the penetration depth A that screening of the mag- 
netic field due to demagnetization is negligible. In a 
cylindrical coordinate system H becomes^: 
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where K° is the kinetic energy in zero external field (ab- 
sorbing the chemical potential) and uj c = eH/m*c is the 
cyclotron frequency. With the typical separation of vari- 
ables, the electron pair for the basis in Eq. consists of 
time-reversed electrons in states (mjk | ) and (—mj—k I) 
with to the azimuthal quantum number in e lmS and k the 
z-wavevector in e %kz . We explicitly write down the BdG 
equations by expanding u(r) and v(r) in terms of the 
eigenf unctions of K° as 
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with L the length of the cylinder. We apply a boundary 
condition 



it(r) = v(r) — with |r| = R, 



(8) 



so that the wavefunction vanishes outside the cylinder. 
This relation only imposes the condition that there are 
no electrons outside the cylinder and does not make any 
assumptions on the coarse-grained superconducting or- 
der parameters as usually treated in Ginzburg-Landau 
theory^. The radial term 
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where J m is the m-th order Bessel function and a m j is its 
j-th zero. The operator K° is diagonal with matrix ele- 
ments K° m]k = h 2 /2m* (a 2 mj /R 2 + k 2 ). The BdG equa- 
tion requires evaluation of matrix elements for (r 2 ) and 



A(r). We will consider only the case of order parameters 
with zero net angular momentum and zero net momen- 
tum along z-axis, namely m + m' = and k + k' = in 
the product u^v^j, for order parameter A(r). A paired 
state with finite net (angular) momentum has higher ki- 
netic energy than the stationary solution, and is therefore 
disfavoredii. With this choice of order parameter, we can 
compute the matrix elements for A(r) and r 2 as 
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The transverse modes (indexed with mj) are decoupled 
from the longitudinal modes (indexed with k) and the 
gap equation is simplified. 

The BdG equations Eq. now become (with Im-.jj' = 



mjk 
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Note that the Zeeman terms |(w+ l)ftiJ c have the same 
sign for u^j and v^j since they represent the amplitudes 
for time-reversed states. The order parameter A(r) is 
sclf-consistcntly expressed by the typical gap equation 
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where the summation range for the eigenstates is over ki- 
netic energies K^jk within a window [— wdi^d] of width 
twice the Debye frequency lud- We use generic parameter 
values suitable for conventional superconductors with T c 
at few Kelvins. A at T = converges at large diameters 
as will be shown and we use the converged value Ao = 3.9 
K as the bulk limit throughout this paper, to* is set to 
the free electron mass, ojd — 100 K and fi = 10000 K. As 
shown below, The Fermi velocity is then vp = 0.55 x 10 8 
cm/sec and the coherence length £o = Hvf /%Ao ~ 350 
nm in the bulk. Results are plotted in dimensionless units 
in this paper. Although we have chosen a particular set 
of parameters, we expect that our conclusion will hold 
qualitatively for conventional low-T c superconductors. 

As the diameter D (= 2R) shrinks, the transverse ki- 
netic energy becomes very sensitive to the boundary 
condition, Eq. (JSJ. Defining a variation of radius 5R, 
SK°/K° = -25R/R. With a small uncertainty m ra- 
dius SR = 1 A around R = 10 nm, for example, 
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5K° ~ 0.02/x ~ 200 K, which is comparable to the Debye 
frequency. We incorporate the effects of variations in the 
wire diameter as a noise in the kinetic energy: 
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where s m jk is a random number uniformly distributed 
in [—1,1] and we set SR < 5 A. Due to the time re- 
versal symmetry of static scattering at the boundary, 
s m jk = s-mj-k- Here, we have partly taken into ac- 
count the radial variations via energy levels, leaving un- 
changed the basis functions and the boundary condition. 
We have sampled equi-spaced z-momenta in the half Bril- 
louine zone [0, 7r] with 200 — 500 points, depending on 
the level of convergence required. Note that experimen- 
tal nanowire samples to date have a significant variation 
in diameter along their length, but this variation is of- 
ten slow on the length-scale of the nanowire width; the 
longest-wavelength variations could be subsumed into an 
adiabatic treatment. 



III. RESULTS 

The density of states of a nanorod can be expressed 
as the sum of one-dimensional densities of states N\ (E) 
displaced by transverse energy eigenvalues E a , 



N(E) = Y,Ni{E-E a ). 



(16) 



Quantized transverse levels strongly affect the density of 
states. Transverse modes are spaced with an average level 
spacing SE = I/N2 with N2 the 2-dimensional density of 
states, N 2 = (irR 2 ){m* /ttH 2 ). Since N% is independent of 
the chemical potential, the qualitative results do not de- 
pend on the particular position of the chemical potential. 
More details will be discussed in subsection IIII Bl How- 
ever, Nx(E — E a ) has a van Hove singularity at E = E a 
and therefore small changes in the chemical potential can 
produce quantitatively different results. Confinement ef- 
fects become strong when the level spacing SE is compa- 
rable to A, i.e., 



R< 



n 2 



*A 



(17) 



This condition becomes -D/£o ^ 0.1 with our parameters 
(and A = 3.8 K). Due to the singularities in N{E), so- 
lutions of the BdG equation are quite sensitive to model 
parameters for D/^q < 0.1, even including a moderate 
smearing from SR. 



A. Large diameters 

Superconductivity under confinement has been well 
studied for samples in the limit of SE <C A and D <C 



A 1,12 - 18 . For cylindrical samples with specular boundary 
conditions, we can explicitly solve Eq. (|14fl in the clean 
limit as a function of external magnetic field along the 
axis, 



ln(A/A ) 



a + \/ a 2 — 1 



■\J \ — a 2 , for a > 1 
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= 0, for a < 1 



where a = h/A, h = evpRH /2c, and Ao is the order 
parameter in bulk without external field. Strassler and 
Wyde»i2 have obtained a result similar to the above equa- 
tion for spherical systems and we follow essentially the 
same derivation for cylindrical systems in Appendix ^ 
In this solution, we make a major assumption by ig- 
noring the spatial variation of the order parameter, i.e. 
A(r) = A. In carrying out an analytic calculation from 
Eq. i|14fl . we make further approximations that the posi- 
tion and momentum in Eq. I|14(l commute and that the 



terms quadratic in the field 



]l m jk in Eq. {T^Jl are 



negligible. When compared with a numerical solution of 
the full BdG equations, the above approximations seem 
reasonable in the large D regime for zero external mag- 
netic field, (until -D/£o ^ 0.57 for our parameter values). 

This clean-limit solution has a field dependence quite 
different from the standard Ginzburg-Landau predic- 
tion. Instead of gradually decaying from the zero-field 
order parameter Ao to zero, A stays constant up to 
H\ = 2cAo/evpR and then drops to zero at the criti- 
cal field H c , (see solid line in Fig.^) 
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exp(3/2)ffi = 0.454 



(19) 



with the flux quantum $0 = hc/2e and the coherence 
length £0 = fovp/irAo. The critical field depends in- 
versely on the diameter, since orbital motions of elec- 
tron in small samples are less influenced by magnetic 
field. H^ in Fig. Q is defined as Hi for R = £0, i.e, 
H^ = 2cA /ev F £ Q . 

We have solved the BdG equations under two different 
conditions: first with the constraint of a spatially uniform 
order parameter and then with the constraint relaxed. 
The quantitative results and their overall line shapes (in 
Fig. ^) are in good agreement with the analytic formula 
Eq. i|TB|) . First comparing the analytic result with the 
uniform-A calculation, we find that the invariant order 
parameter up to H = Hi is well reproduced. The numer- 
ical A deviates downward from the analytical formula at 
higher values of H, resulting in smaller critical field H c . 
We attribute the discrepancy to underestimation of the 
external field in the analytic solution due to ignoring the 
A 2 term. Further discussion is given at the end of this 
subsection. 

An important qualitative deviation of the full BdG so- 
lution from the uniform-A results comes at high fields 
near H c . The critical field for the full BdG solution (solid 



circles in Fig. ^) extends to high fields, with A decaying 
much slower than the usual square-root dropiS of the 
uniform-A solution. It is clear that this phenomenon 
must originate from the transverse spatial variation of 
the order parameter, which we will discuss below in more 
detail. 

Although the average order parameter remains con- 
stant for < H < Hi (left to the dashed line in Fig.QJ, 
the density of states reveals a closing of the excitation 
gap at the field Hi (dashed curve in Fig. EJ, where the 
spectral weight at the chemical potential becomes finite. 
The density of states in the gapped region continues to 
increase until H = H c . In the clean limit, the elec- 
tron trajectories are not perturbed by impurity scattering 
and their angular momentum L couples to the external 
field as — (e/2m*c)L ■ H. This coupling gives the term 
— ^mhw c in Eq. (|13|l and contributes to the quasi-particle 
excitation energy -E*^, 

E k m] = V / ^ 2 + A2 + I(m + l)*u> c . (20) 

These angular momentum contributions dominate over 
spin contributions and the excitation spectrum goes gap- 
less when A ~ (e/2m*c)(L ■ H) ~ evpRH/2c. In con- 
trast, when impurity scattering dominates, the electron 
trajectories are disrupted and the quasi-particle energy is 
no longer in the form of Eq. I|2l)fl . In this limit, the quasi- 
particle energy suffers significant level broadening and 
experiences shifts only from the spin Zeeman termsSi. 
Therefore, away from the clean limit the density of states 
retains the form of a conventional gapped superconduc- 

tojis. 

At zero field, the order parameter is nearly constant 
over the cylinder (see Fig. except for small oscilla- 
tions and a Gibb's phenomenon at r = R. The rapid os- 
cillations have a wavelength proportional to 1/vf and an 
amplitude that diminishes for larger diameters. There- 
fore we expect that these oscillations will be averaged out 
on larger length scales; they are not important for large- 
diameter systems. In such regime, we correctly reproduce 
the boundary condition commonly used in the literature, 
namely a vanishing normal derivatives of the supercon- 
ducting order parameter at the surface. The order pa- 
rameter does not change until H = H\. As H exceeds 
Hi, some quasi-particle energies are pushed below zero 
in Eq. 12U|) and contributions from these excited quasi- 
particle states reduce the order parameter in Eq. (fH|) by 
changing the signs of their contributions in the statistics 
factor, 1 — 2/(£* y ). In addition to an overall reduction 
of A(r) under external field, A(r) also changes slope, 
with a distinct knee that moves towards r = with in- 
creasing field (see arrows in Fig. At a radius r, the 
angular momentum is of order m*rvF and the order pa- 
rameter begins to be suppressed when A ~ m*vprH/2c 
or r ~ 2cA/evpH . This simple argument should be 
taken with caution, since the semi-classical approxima- 
tion of treating position and momentum as commutable 
in Eq. I|14|l becomes worse when there is a strong spa- 
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FIG. 1: Averaged superconducting order parameter in a 
cylindrical sample of large diameter (D = 200 nm, -D/£o = 
0.57) with axial external magnetic fields. Normalized order 
parameter to the zero field value A/Ao is plotted as a function 
of dimensionless field HR/(H^o) with = 2cAo/(evF^o)- 
Overall agreement of the analytic formula Eq. 1181 and the 
numerical results are good. Solution of Bogoliubov-de Gennes 
(BdG) equation in the clean limit remains constant until 
Hi = 2cAo/cvfR- Full solution of the BdG equation (solid 
circles) has larger order parameters near the critical field than 
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FIG. 2: Density of states as a function of external field. 
Superconductor becomes gapless at H — Hi (dashed line) in 
the clean limit well before the critical field H c , due to the 
coupling of orbital angular momenta to external field. The 
external fields in the plot are from HR/(H^o) = (thick 
line) to 1.83 with equal intervals between the curves. 

tial variational At larger radius r, the energy difference 
between the angular momenta (with mh ~ m*rvF) in 
an electron pair exceeds the pairing energy and therefore 
the pair becomes depaired. 

We emphasize that the pronounced radial dependence 
of order parameter is related to the coupling of orbital 
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FIG. 3: Spatial variation of order parameter A in the large 
diameter limit (-D/fo = 0.57). For external fields of < H < 
Hi, A remains flat. As the field increases from Hi, A drops 
with a knee which progresses towards r = 0. 
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FIG. 4: The order parameter A as a function of the radius 
at small diameter D/£ = 0.057 (SR = 2 A). The spatial 
variation of A is much stronger than for the large diameter 
wires of Fig. [3] 



B. Small diameters 



angular momentum to the external field, rather than the 
term ^m*u 2 r 2 in Eq. JSJ, which makes the superconduct- 
ing order parameter more massive. To compare these two 
contributions, we consider the Ginzburg-Landau theory 
with an order parameter of zero total angular momen- 
tum, as is usual in the literature^; (see appendix lB|l . 
The order parameter couples to the external field only 
through ^m*uj 2 r 2 and does not have any information 
about the angular momenta of constituent electrons in 
Cooper pairs. The spatial dependence of A(r) arising 
from the term (2e / he) 2 A 2 turns out to be much weaker 
than the angular momentum coupling in the BdG treat- 
ment. The resulting critical field Eq. l|B5|) is almost iden- 
tical to the angular momentum depairing result Eq. (|19H • 
as both are derived in detail in the appendix. 

It is interesting to note that, although the orbital cou- 
pling and the A 2 -coupling produce very close critical 
fields (see Eqs. (jl9IJj5|l ) in the absence of the other, they 
do not act additively when both are present. When the 
A 2 -teim is included on top of the orbital coupling, its 
effect becomes significantly smaller than in Eq. I|B5|) be- 



cause the orbital effect produces fast-decaying |A(r)| 2 
near H = H c . For example, the curve labeled with 
HR/(H^ ) = 2.16 in Fig. El has \A(r = R)\ 2 /\A(r = 
0)| 2 as 1/7. Since the A 2 ~ H 2 r 2 term couples most 
strongly at large r, the actual A(r) 2 |V'( r )| 2 coupling is 
much smaller than the Ginzburg-Landau theory. For 
instance, a reduction of the effective A 2 -coupling by 
half results in the reduction of the order parameter by 
+ 0.5) as 82% in terms of the Ginzburg-Landau 
theory when the both couplings are naively added, which 
nearly matches the discrepancy in Fig. ^ between the 
Eq. i|18|) (thin line) and the constant-A (open circles) 
results. 



When the diameter shrinks sufficiently that the trans- 
verse level spacings SE exceed the order parameter A, the 
density of states on the energy scale of A becomes spiky 
and spatial structure arising from the transverse modes 
begins to show up in the radial dependence of A(r). Since 
the level spacing of the transverse modes is inversely pro- 
portional to the effective mass (SE = h 2 /m* R 2 ), the ef- 
fect of the discrete levels will be strong for systems of 
small effective mass or low transition temperature. Fig.0] 
shows the spatial variation of the order parameter in this 
regime. A close examination reveals two characteristic 
length scales in A(r). The shorter length scale (with 
rapid oscillations more apparent at large radius) is given 
in terms of the Fermi velocity, i.e., Ar ~ 2 r Kh/m*VF- As 
the Fermi velocity (or the carrier density) grows, A(r) 
oscillates more rapidly. 

Apart from the structures corresponding to the Fermi 
wavelength, there are more interesting and slowly varying 
spatial modulations, particularly near r = 0. Although 
these modulations also appear for large diameters (see 
Fig. |3J) , their relative importance grows in the smaller di- 
ameter wires. This spatial structure arises from the small 
number of transverse modes within the energy window 
[-Wrj,Wfl]. As can be seen from Eqs. ((141 IA3|) . the en- 
ergy levels close to the Fermi energy contribute strongly 
to the order parameter A. Since the density of states is 
peaked at the transverse energy levels (see Eq. I16f) . the 
resulting order parameter has larger amplitude for the 
states with h 2 a 2 n j/2m*R 2 as /z, k z as and displays the 
spatial characteristics of those transverse modes. As il- 
lustrated in Fig. 03 the states within the fc-space shell of 
area 2nkpSk± ~ 47rm*urj contribute most strongly to 
the order parameter. Although we have used sharp en- 
ergy cut-offs at lod, they are not expected to impose a sig- 
nificant quantitative change because the weight factor in 
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FIG. 5: Schematic 2-dimensional phase space with zero ; 
ial wavevector {k z = 0). Transverse modes within the si 
of thickness Sk± contribute strongly to the order paramet 
The area of the shell 2nkp5k± ~ 4nm*ujD does not depe 
on the choice of Fermi energy. 



the gap equation, u a v a ps A/2y / e^ t + A 2 (cf. Eq. 

is small near the cutoff. The ratio of u a v a for e a = to 

e a = is approximately wrj/A w 25. 

The shape of the order parameter is determined by 
which states happen to fall into the fc-space shell. For 
instance, if states of (m = 0) are absent in the shell, then 
the amplitude A(r ~ 0) is depleted, since only the Bessel 
functions of Joj have non- vanishing values at r = 0. We 
caution that the order parameter is not necessarily zero 
at r = 0, since there are also states (mjk) with finite k. 
As the diameter decreases, the transverse states become 
more sparse in the energy shell and the spatial structure 
becomes more pronounced. In contrast, large-diameter 
wires have many contributing {mj, k = 0) states and the 
spatial variation averages out. 

Averaged order parameters A at zero magnetic field 
are plotted as a function of diameter in Fig. H3 The tem- 
perature is fixed at 0.2 K and 3 different radial smearings 
SR = 0, 2, 5 A are used. The filled circles are solutions 
for full BdG equations and the open circles impose the 
constraint of a constant A. Regardless of SR and the 
spatial constraint, A/A converges to 1 at large diame- 
ter. A fluctuates considerably as the diameter decreases, 
with more scatter for smaller radius smearing. These 
variations arise from the sharp van Hove singularities in 
the density of states at transverse eigenvalues. Interest- 
ingly, the full BdG order parameter is consistently larger 
than the constant A solution. The gap equation, Eq. I|14|) 
or (|A3|) , becomes particularly simple for a constant A at 
T = 0, 



1 

V 



Vuj 2 + A 2 



, with p D (u) = y^£(e a - 



( 21 ) 

where e a is the non-interacting eigenvalue (absorbing 
chemical potential). Although the density of states 
Pd(&) at diameter D has fluctuations due to discrete 
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FIG. 6: Average order parameters A as a function of diameter 
for different radius smearings SR. Filled circles are the full 
BdG solution while empty circles represent the constrained 
case where spatial variation of A is disallowed. Confinement 
effects appear at D/£o ~ 0.1, i.e., when the transverse energy 
level spacing 5E w A. A converges to about 4 K as the 
diameter increases, regardless of SR. In the full solution, A 
increases slightly as D/£o drops below about 0.1. Compared 
to the constant-A behavior, this enhanced order parameter 
takes advantage of the spatial variation in the BdG solution. 



transverse energy levels, it averages to the bulk limit. 
Therefore, statistically, A is expected to fluctuate about 
the bulk value An . Particular lineshape of the curves (of 
empty circles) in Fig. [S] is due to limited diameter sam- 
pling and the fine-tuning of model parameters. As the 
calculation of SR = 0.5 nm indicates, stronger broad- 
ening suppresses the fluctuation of the order parameter 
which will converge to An down to small D. When the 
condition of the uniform A is relaxed, the A has the free- 
dom to peak in regions with a higher density of electronic 
states, thereby increasing the condensation energy. For 
all SR shown here, the enhanced order parameter is most 
evident for D/£q smaller than about 0.1, the diameter 
regime where the transverse level spacing becomes com- 
parable to A. For D/£q of 0.05-0.1, the enhancement is 
roughly 10 — 20 %. This is consistent with larger T c in 
small samples, as is often observed in thin filmsSiSii 2 ^. 
This tendency has been attributed to softening of surface 
phonons in samples of small dimension. The spatial vari- 
ation of the order parameter A from confinement effects 
may also contribute to this trend. 

The order parameter A versus magnetic field plotted 
in Fig. shows shoulders that also reflect the discrete 
nature of the transverse modes. The overall shape of 
the curves is similar to that for large diameters (see 
Fig. ^) . A remains constant until the depairing field 




FIG. 7: The order parameter A as a function of the external 
field at a small diameter D/Co = 0.057 - 0.086 (D = 20 - 30 
nm, SR = 2 A). A displays several shoulders as the field in- 
creases and then vanishes abruptly, unlike the large diameter 
case of Fig. This structure reflects the discrete nature of 
the transverse modes. The sudden cut-off of A at the critical 
field is due to an absence of small angular momentum trans- 
verse states near the chemical potential. The dashed line is 
Eq. iflSjl. 

in Eq. l)A3j) becomes comparable to A. As the field in- 
creases further, distinctive shoulders appear. As shown 
in Eq. (|A3|) . until the condition \J e 2 a + A 2 < eha is sat- 
isfied for any non-interacting state a, the thermal factor 
1 — 2/ does not change, so the gap equation yields the 
same A. At large di ameters, the number of depaired 
states with ^/e 2 + A 2 > e/ lQ increases gradually and A 
therefore varies steadily with the external field H. How- 
ever, at small diameters, the density of states has peaked 
structures with van Hove singularities separated by an 
energy spacing of h 2 /m*R 2 . Since only a few transverse 
energy peaks are available in the energy window, the van 
Hove singularities have a stronger influence on the order 
parameter. Therefore kinks begin to appear in the field 
dependence of A. As D grows, these discrete structures 
smooth out, as in Fig.^ 

The independence of the critical field H c is shown in 
Fig. El Electron orbits in more tightly confined spaces are 
less influenced by magnetic field, because of the smaller 
depairing contribution (e/Hc) J d£ ■ A in the phase of 
wavefunction. The filled circles are BdG solutions with 
5R = 2 A and the dashed line is the semi-classical (i.e. 
l/R) solution of the BCS equation Eq. (fl^|l . The inset 
magnifies the small diameter regime. The BdG solution 
follows the l/R trends well down to about £>/£o ~ 0.02. 
As D gets smaller, H c fluctuates substantially, but fol- 
lows the overall l/R behavior surprisingly well. One of 
our finding is that the prediction of the Ginzburg-Landau 
theory agrees remarkably well with the microscopic so- 
lution of the BCS equation, down to small diameters 
(D <C Co)- This conclusion may change with an inclusion 
of strong interaction effects, such as increasing Coulomb 



FIG. 8: The critical field vs rod diameter. The critical field 
H c varies inversely proportional to the diameter D. The so- 
lutions of the full BdG equation follow the general trend of 
Eq. lO, which was obtained from the BCS equation under 
the constraint of a spatially constant order parameter. The 
full BdG solution satisfies the 1/7? prediction down to fairly 
small diameters, far smaller than coherence length, although 
they fluctuate strongly about the relation Eq. (11911 . 

interaction in strong confinement at small diameters. 

Finally, we mention that temperature dependence of 
the order parameter does not show significant deviation 
from the BCS results^. The relation A/k B T c = 1.764 
holds to high accuracy for a wide range of diameters at 
zero external field. 



IV. CONCLUSIONS 

We have studied the dynamics of the transverse de- 
grees of freedom in superconducting nanowires. As the 
confinement dimension shrinks and the level spacing be- 
comes comparable to the order parameter A, the dis- 
crete nature of the transverse modes shows up in a spa- 
tial variation of A. In the clean limit, electronic angular 
momenta are conserved and couple strongly to magnetic 
field to shift the quasi-particle energy levels. This ef- 
fect shows up as distinct shoulders in the response to 
an external field. In a confinement scale comparable to 
or larger than the superconducting coherence length, su- 
perconducting wavefunctions satisfy the usual boundary 
condition for normal derivative dip / dn — for supercon- 
ducting order parameter ip. Under zero magnetic field, 
■0 remains constant throughout the sample, except for 
small and rapid oscillations. With finite external field, ip 
adjust to the vector potential with much stronger spatial 
variations than predicted in Ginzburg-Landau theory. 

Our results are relevant to clean-limit samples with 
inclusion of level broadening effects introduced by uncer- 
tainties in diameter. It is useful to compare the results 
with the dirty-limit theories£2*i2*iii. Although detailed 
mechanisms for both limits are different, both systems 
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display gapless behavior. The critical fields for the dis- 
appearance of order parameters, H c , behave quite differ- 
ently for the two limits. In the dirty Urnffci£(£/£o <c 1), 
H c becomes very large (before the spin-Zeeman depair- 
ing effect dominates^), with H C {1) ~ H c ^£q/£ > H c 
with H c given in Eq. I|19|l for the clean limit (£/^q 
1). For the critical field Hi where the excitation spec- 
trum first becomes gapless, Hi/ H c \i = o = 0.954 while 
Hi/ H c \z =00 — 0.389. Therefore, the clean limit remains 
gapless for a wide part of the magnetic field range com- 
pared to only 4.6% of the dirty limit. While the on-set of 
the order parameter suppression and the disappearance 
of the excitation gap happen simultaneously in the clean 
limit, the closing of the excitation gap in the dirty limit 
happens only when the superconducting order is already 
suppressed significantly. 

Electrodeposition into nanoporous membranes such 
as polycarbonate or anodic alumina^ can yield single- 
crystal metallic nanowires from several different super- 
conducting metals (tin, lead, etc.). Such systems may be 
able to access the clean limit in which the phase infor- 
mation of definite angular momentum states is conserved 
and orbital-derived level-shifts under magnetic field be- 
come substantial. Since the effects of discrete levels begin 
to appear when the level spacing becomes comparable to 
the order parameter A ( cf. Eq. (|f 7[0 . systems of smaller 
A will exhibit stronger confinement effects at a given wire 
diameter. Systems with small band mass m* will have 
similarly strong confinement effects. 

The nature of confined superconductivity in single- 
crystalline metallic nanorods could perhaps be verified 
most clearly by the gapless spectrum that appears at 
magnetic fields smaller than the critical field. The small 
quasi-particle excitation energy here could result in very 
interesting physics, e.g. in specific heat measurements 
at temperatures below T c under external field. Gapless 
superconductivity under an external field could also en- 
hance phase-slip rates as reflected in the electrical re- 
sistivity. Due to the low quasi-particle excitation ener- 
gies, thermal or quantum fluctuations will overcome the 
condensation energy more easily. Gapless excitations for 
normal electrons could contribute to a finite residual re- 
sistance in the presence of strong scattering at the con- 
fining surface. 



APPENDIX A: SEMI-CLASSICAL SOLUTION 
FOR ORDER PARAMETER A IN CYLINDER 

We derive the analytic solution to the gap equation 
Eq. I|f 4f> with the constraint of a uniform order parame- 
ter A. The relevant regime is when the cylinder is large 
enough that the transverse level spacing is much smaller 
than A and the external field is not too close to the criti- 
cal value. Similar results have been derived for spherical 
systems using a Green function technique^ With the as- 
sumption of a constant A, one approximates [r, p] w 0, 
treating r and p as independent variables. We ignore the 
Zeeman term from spins, /i^er ■ H in Eq. J3J), since it is 
negligible in the regime of interest. The term quadratic in 
the field, e 2 A 2 /2m*c 2 , is also dropped, as discussed be- 
low. Under these conditions, the non- interacting Hamil- 
tonian becomes Hq = £k — eh, where £k = \m*v 2 — /i 
and 



-L H. 



.. (Al) 

2m* c v ; 

Solving the BdG equation Eq. Q for a constant A, one 
obtains 

E u = Jel + A 2 - e h , and M k w k = - . ,f ■ (A2) 

2 v £ k + A 

With these approximations, u^v^ is independent of the 
field H. The gap equation can be written as 



1 = V 



E 

k 



l-2/(y^+¥-e fe ) 
2,/4TA^ 



(A3) 



The effect of the magnetic field is reflected only in the 
Fermi-Dirac function /. To evaluate the angular mo- 
mentum summation, we write e/,. = (e/2m*c)L • H = 
(e/2c)(r x v) • H = (e/2c)(H x r) • v. Now, with 
a fixed r and H and an isotropic distribution of v, 
e/j = (evF /2c) rH cos <p with the angle ip = Z(H x r,v) 
and |v| w vf- The summation on the right-hand side of 
the gap equation Eq. (|A3|I becomes 



deN 



2irdr \ 



dQ v 1 - 2/ 



irR 2 J Air 2%/e 2 + A 2 ' 



(A4) 



where N is the density of states and fi„ the solid angle 
for v. If e 2 A 2 /2m*c 2 is much smaller than the integral 
limit u>d, it can be absorbed in the chemical potential 
with little change in the integral J de. The gap equation 
becomes 



f = N V 



de duu 
Jo 



dft 



I - 2/(Ve 2 + A 2 - hufi) 
Ve 2 + A 2 
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(A5) 

with u = r / R, \i = cos ip and h = (evp /2c)RH . For 
H = 0, we have 
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f = N V 



= N Vln 




2/( v /^TA 2 ) 



A 2 , 



(A6) 
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Since / = for \J e 2 + A 2 > hu/j, in the T = limit, we 
can perform the integral in Eq. (|A5I) 



de I duu 
Jo 



1 /(Ve 2 + A 2 - feu/*) 



1 



1 

2^ 



111 



Ve 2 + A 2 
3 



1 



y/T 



fora = )i/A> 1, and for a < 1. For a > 1, Eq. 
becomes 



N V 




Substracting Eq. i|A6(l from Eq. (|A7|) and using > 
Ao, we obtain Eq. JTHJ. The order parameter A be- 
comes gapless as one or both of the states in the time- 
reversed pairs in Eq. (|A2fl are pushed out of the Debye 
frequency. The external field influences T c by the statis- 
tical factor, not directly through the field dependency in 
the non- interacting density of states, as can be inferred 
from Eq. (|A5|) . The density of states over the interaction 
window of the Debyc frequency remains nearly the same 
due to a balance between the outflux and the influx of 
non-interacting energy levels shifted by the external field. 



H c is the bulk critical field, ip — tp f (|?Ao| 2 = ol/(3), 
lu c = eH/moc, mo is free electron mass, and £ is the co- 
herence length. The Coulomb gauge is used for the vector 
potential A = ^rH8. For radius R much smaller than 
the penetration depth, the screening of external fields is 
negligible, as assumed in the main text. With constant /, 
we obtain the average free energy after integration over 
the cylindrical nanowire, 



Hi 
8tt 



1 / m uj c RC 

2 V » , 



(B3) 



We have also numerically minimized the free energy with 
a spatially varying /, but the spatial variation was far 
less significant (about 1 percent of the average order pa- 
rameter for R/£ = 0.5) than that of the Bogoliubov-de 
Gennes results discussed in the text. Therefore, the fol- 
lowing approximate value for the critical field is quite 
accurate. Minimizing Eq. I|B3|> over / gives 



/ = 



1 / m uj c R(, 
2\ h 



"I 1/2 



(B4) 



APPENDIX B: GINZBURG-LANDAU THEORY 

Here we consider only a stationary superconducting 
order parameter [L z — 0, k z = 0) where the Ginzburg- 
Landau functional for the order parameter ip reads 



F = - a |V,| 2 + ^| 4 



^ .2e A . 

V - i— A ip 

he 



(Bl) 



-Iff 



v- 



imQLU c rb 



For vanishing /, the critical field is 



He = —t|= 0.450-| , 

7T t;R t;R 



(B5) 



.(B2) with $ = hc/2e. 
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